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ABSTRACT 

The  practice  of  forming  electrical  conduction  paths  in  an  insulating  material  by  filling 
cylindrical  holes  with  molten  metal  can  result  in  high  residual  stresses  when  the  metd 
cools.  Residual  stress  is  greatest  near  the  metal-insulator  interface,  and  stress  relaxation 
by  mecins  of  de-adhesion  is  possible.  Another  failure  mode  that  poses  greater  practical 
difficulties  is  the  growth  of  cracks  along  paths  which  spiral  away  from  the  interface  into 
the  brittle  materid.  Such  cracks  may  occur  singly  or  in  pdrs,  and  their  lengths  can 
be  sufficiently  great  to  provide  links  with  adjacent  conduction  paths.  Such  cracks  are 
considered  from  the  fracture  mechanics  point  of  view.  The  residud  stress  field  is  relaxed  by 
the  growth  of  spird  cracks  which  are  modeled  as  continuous  distributions  of  dislocations. 
It  is  assumed  that  these  cracks  grow  so  that  the  stress  state  on  the  prospective  fracture 
plane  just  ahead  of  the  crack  tip  is  purely  tensile.  The  paths  are  determined  by  means  of 
an  incrementd  numericd  procedure. 


INTRODUCTION 

The  question  of  reliability  of  integrated  circuits  is  concerned  with  the  predictable 
sustained  performance  of  some  electronic  function.  Even  though  the  capacity  for  crirrying 
mechanical  loads  may  not  be  relevant  to  this  function,  fdlures  of  devices  are  commonly 
due  to  mechanical  effects.  These  failures  are  driven  by  residud  stresses  arising  from  fab¬ 
rication  processes,  often  carried  out  at  elevated  temperature  and  involving  combinations 
of  materials  with  very  different  thermd  tmd  mechanical  properties.  The  fracture  problem 
studied  here  is  of  this  kind. 

In  fabricating  an  integrated  circuit  with  a  complex  planar  or  multi-planar  structure, 
it  is  sometimes  necesstuy  to  pass  conduction  lines  directly  through  the  brittle  substrate 
of  a  device.  This  may  be  done  by  boring  a  hole  through  the  thickness  of  the  substrate, 
which  is  assumed  to  be  inilidly  stress  free.  The  hole  is  then  filled  with  molten  metd,  say 
duminum  or  copper,  and  the  molten  metd  adheres  to  the  surface  of  the  hole  as  it  cools 
and  solidifies.  As  a  result,  an  elastic  mismatch  extensiond  strain  Ae  can  arise  between 
the  metd  inclusion  in  the  cylindried  hole  and  the  surrounding  brittle  materid. 

If  the  brittle  substrate  is  thick  compared  to  the  radius  of  the  hole  (see  Figure  1), 
then  the  mechanied  fields  over  most  of  the  length  of  the  inclusion  are  essentidly  two 
dimensiond  plane  strdn  fields.  This  is  the  case  studied  here.  The  tensile  stress  in  the 
residual  stress  field  tends  to  be  greatest  at  tiie  bimaterid  interface.  Indeed,  fractures 
typicdly  nucleate  at  the  interface.  However,  it  is  often  observed  that  fractures  grow  away 
from  the  interface,  and  they  tend  to  spird  Eiround  the  cylindrical  inclusion,  as  suggested 
in  Figure  1.  Furthermore,  the  detailed  shap>e  of  the  spird  appears  to  be  quite  sensitive 
to  the  relative  properties  of  the  two  materids  involved.  The  purpose  here  is  to  examine 
the  spird  cracking  process  from  the  perspective  of  brittle  fracture  mechanics. 

'I'he  stress  field  of  a  two  dimensiond  elastic  solid  contdning  a  sharp  crack  is  nom¬ 
inally  square  root  singular  at  the  crack  tips.  The  amplitude  of  this  singularity  is  cdled 


91-00871 

[i  iiiii  ttiti  .-T:  "  ■ 


Presented  at  the  1991  MRS  Spring  Meeting,  Anaheim,  CA  in  Symposium  H: 
Mechanical  Behavior  of  Materials  and  Structures  in  Microelectronics 


91  &  SO  030 


Fig.  1.  Sdicinatic  diagram  of  cyliii- 


the  stress  intensity  factor  and  the  spatial  dependence  of  the  asymptotic  field  is  univer¬ 
sal.  The  tendency  for  crack  growth  is  commonly  expressed,  as  a  postulate  on  material 
behavior,  in  terms  of  critical  or  material-specific  values  of  the  stress  intensity  factor.  Of 
particular  concern  here  is  the  connection  between  the  direction  of  crack  growth  and  the 
stress  intensity  factor.  For  plane  strain  deformation,  the  crack  opening  near  the  tip  cem 
be  symmetric  with  respect  to  the  local  plane  of  the  crack  (mode  I),  antisymmetric  (mode 
II),  or  some  combination  of  the  two  (mixed  mode).  It  has  frequently  been  observed  that 
cracks  in  brittle  materi2ils  tend  to  advance  in  such  a  way  as  to  maintain  a  symmetric 
deformation  field  (mode  I)  near  the  moving  crack  tip  [1,2).  This  additional  postulate 
on  fracture  behavior,  sometimes  called  the  principle  of  local  symmetry,  is  adopted  as  a 
criterion  to  predict  the  path  of  crack  advance  around  the  cylindrical  inclusion. 

The  analysis  proceeds  in  several  steps.  First,  the  residual  stress  field  due  to  the 
thermal  mismatch  strain  is  estimated  in  the  absence  of  a  crack.  Then,  a  crack  is  assumed 
to  initiate  at  the  interface  and  to  grow  2iround  the  inclusion.  A  numerical  algorithm  is 
developed  foi  extracting  the  values  of  the  mode  I  and  mode  II  stress  intensity  factors, 
denoted  by  K /  and  Kn,  zs  the  crack  gets  longer.  The  same  procedure  h2is  been  applied 
for  straight  cracks  in  a  homogeneous  body  in  [3].  The  algorithm  is  sufficiently  general  to 
analyze  crack  growth  along  a  path  which  is  not  specified  in  advance  but  which  is  selected 
incrementally  by  the  local  symmetry  condition  Ku  =  0.  Results  are  presented  for  cases 
when  the  stiffness  of  the  inclusion  is  equal  to,  greater  than  and  less  than  the  stiffness 
of  the  surrounding  brittle  elastic  material.  Similar  calculations  restricted  to  the  case  of 
homogeneous  bodies  have  been  presented  in  (4)  and  [5]. 


INITIAL  RESIDUAL  STRESS 


Suppose  that  a  cylindrical  hole  of  radius  a  exists  in  a  stress  free  block  of  brittle 
material  which  has  dimensions  large  compared  to  a.  The  hole  is  filled  with  molten  metal 
which  adheres  to  the  surface  of  the  hole  as  it  cools.  Suppose  that  an  extensional  mismatch 
strain  Ae  results  from  the  temperature  change  emd/or  the  difference  in  thermal  expan¬ 
sion  coefficients.  The  axially  symmetric  stress  distribution  in  the  surrounding  material, 
referred  to  polar  coordinates  centered  on  the  hole,  is 

0 _ 0  _  _  fiC  o  _  „ 

p„(l +  !/„)  + p(l-«/„)r2  -  2x(l-,.)r»’ 

where  p  and  i/  are  the  shear  modulus  and  Poisson’s  ratio,  respectively,  of  the  brittle 
insulator  and  pm  and  Vm  are  the  corresponding  constants  for  the  inclusion.  The  stress 
(1)  is  to  be  relieved  through  the  process  of  spiral  cracking.  All  computed  results  for  stress 
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introduced  in  (1).  When  the  elastic  properties  of  the  inclusion  and  the  surrounding  brittle 
material  are  approximately  equal  then  C  has  a  value  about  ten  times  the  elastic  mismatch 
strain  Ac.  All  values  of  stress  intensity  factor  calculated  here  are  normalized  by 


A" 


fiC\/a 

2(1 


(3) 


To  give  an  indication  of  the  value  of  this  normalizing  parcuneter  for  a  particular  sys¬ 
tem,  consider  an  inclusion  with  radius  a  ~  0.1mm  and  mismatch  strain  of  0.001  in 
a  material  with  shear  modulus  of  about  30  GPa.  These  numerical  values  imply  that 
K'  ~  2MPav/^. 


THE  SUPERPOSITION  SCHEME 

Consider  any  curve  in  the  plane  that  is  a  prosp>ective  path  for  a  spiral  crack,  and 
denote  this  curve  by  £.  A  pair  of  orthogonal  unit  vectors  e»,  Cn  aligned  with  the  tangent 
and  normal  to  C  are  introduced,  as  shown  in  Figure  2.  The  variable  s  is  the  arclength 
along  C  measured  in  the  direction  of  Ci  from  some  arbitrary  point.  The  angle  ^(s)  is  the 
angle  measured  counterclockwise  from  the  z-axis  to  the  tangent  line  of  C  at  s.  For  any 
state  of  stress  in  the  body  expressed  with  respect  to  the  polar  coordinate  directions  or 
the  underlying  z,y  rectangular  coordinate  directions,  the  normal  traction  t„(s)  and  the 
shear  traction  T((s)  associated  with  the  normal  vector  en  at  any  point  eilong  C  can  be 
expressed  as  the  complex  quantity 

r„(s)  +  ir,(s)  =  i  +  cr,,)  +  [1  (cr,»  -  <T**)  + 

(4) 

=  i  {<7rr  +  +  [5  {(Ttt  -  ‘^rr)  +  iOrs] 

where  t  =  %/— 1.  In  this  expression,  it  is  understood  that  the  coordinate  s  on  £,  coordi¬ 
nates  z,y  in  the  plane,  and  coordinates  t,6  in  the  plane  represent  the  same  point. 

In  terms  of  the  representation  (4)  and  the  residual  stress  (1),  the  traction  to  be 
relieved  by  crack  growth  is  <t®(s)  +  ifff(a).  The  process  of  crack  growth  is  simulated  by 
introducing  a  continuously  distributed  array  of  Volterra  elastic  dislocations  for  which  the 


s 


distribution  is  selected  to  cancel  the  initial  residual  stress.  This  is  a  common  device  for 
analyzing  fracture  phenomena  in  brittle  materials  [6].  The  fundamental  building  block 
in  the  procedure  is  the  expression  for  the  complete  stress  field  due  to  a  single,  isolated 
dislocation  with  arbitrary  Burgers  vector  at  an  arbitrary  point  in  the  elastic  body  being 
analyzed.  In  the  present  case,  the  body  is  the  whole  plane  with  elastic  constants  E  and  v 
containing  an  embedded  circular  inclusion  of  radius  a  with  elastic  constants  Em  and  Um- 
This  stress  distribution  can  be  found  without  approximation  as  outlined  in  Appendix 
I.  Features  of  particular  significance  for  this  development  are  that  each  stress  component 
is  linear  in  the  components  of  the  Burgers  vector  of  the  dislocation  and  that  each  stress 
distribution  has  a  pole-type  singularity  at  the  point  coinciding  with  the  dislocation  in 
the  plane.  Consider  once  again  the  arbitrary  curve  C.  Supp>ose  the  dislocation  is  located 
at  a  point  on  that  curve  identified  by  a  particular  value  of  arclength,  say  s,  and  that 
the  Burgers  vector  of  that  dislocation  is  denoted  by  components  hf{a)  and  6y(s).  Now 
consider  any  other  point  on  C  located  by  means  of  a  value  of  arclength  q.  The  features 
of  the  stress  distribution  noted  above  imply  that  the  traction  components  (r„  and  <T|  at  q 
due  to  the  dislocation  at  a  can  be  written  in  the  form 


0'n(9i^) —  [^nr(9»  d"  ^)^y(^)] 

3  —  q 

=  — ^[E,x(9,s)6,(s)  +  E,y(9,s)fcy(s)] 
.s  ^ 


(5) 


The  steps  required  to  evejuate  the  four  functions  E„z(q,s), . . .  appearing  in  (5)  for  any 
values  of  s  and  q  are  given  in  detail  in  Appendix  II.  These  functions  play  a  central  role 
in  the  solution  procedure. 

The  stress  fields  discussed  above  zire  now  generalized  from  the  case  of  a  single,  isolated 
dislocation  with  finite  Burgers  vector  to  the  case  of  a  continuous  distribution  of  disloca¬ 
tion  with  infinitesimal  Burgers  vectors.  This  transition  is  effected  simply  by  making  the 
replacements 

bz{^)  — *  ^z{s)ds  and  fcy(s)  — *  ipy{s)ds  (6) 

and  then  allowing  the  arclength  s  to  range  over  values  s_  <  s  <  s+  on  £.  With  this 
interpretation,  i/>i(s)ds  represents  the  net  Burgers  displacement  in  the  t-direction  of  a 
dislocation  smeared  out  over  the  infinitesimal  interval  from  s  to  s  -f  da.  In  other  words, 
0i(s)  is  the  gradient  along  £  of  the  discontinuity  in  the  i-component  of  displacement 
across  £  at  s.  If  the  discontinuity  in  displacement  Ui  at  a  is  denoted  by  [uij(s),  then 


[“rl(9)=j^  xl^r(s)ds,  [Uy]{q)  =  J^  ds 


(7) 


where  it  is  tacitly  assumed  that  the  crack  is  closed  without  net  dislocation  at  s  =  s_.  If 
the  same  is  true  at  s  =  s+  then  if-z  emd  i/>y  are  subject  to  the  closure  conditions 


J  xl>z{s)ds  =  0,  J  ipy{s)ds  =  0 


(8) 


The  functions  ipz  and  xpy  are  expected  to  be  square-root  singular  at  the  crack  ends  [6]  and 
to  be  asymptotically  proportional  to  the  elastic  stress  intensity  factors  there.  In  terms 
of  the  mode  I  emd  mode  II  stress  intensity  factors  discussed  in  the  Introduction,  this 
asymptotic  relationship  is 


Kn  -1-  iKj  =  ^hm^  v/2jr(3+  -  a)  [V>.(s)  +  »V>,(s)]e~''*^*> 


(9) 
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For  the  continuously  distributed  dislocations,  the  traction  at  axiy  point  9  on  £  is 
still  given  by  (5),  provi«led  that  the  expressions  on  the  right  side  of  (5)  are  summed  over 
all  values  of  s.  This  is  accomplished  by  m^lking  the  replacements  (6)  in  (5),  and  then 
integrating  over  the  full  range  of  s  on  £.  The  result  is 


/■*+  fig 

<^n(9)=  /  (E  nx  (9,  •s)t^’r(«)  +  E„,(9,  s)V’»(3)]  - 

Ja-  5  — 

=  /  [Sii(9, 3)il>r(s)  +  E,y(9,  s)^y{s)] 

Ja^  5  —  ^ 


(10) 


Thus,  if  the  curve  representing  the  crack  line  £  is  known,  and  the  distribution  of  crack 
opening  t/’i(3),  along  the  curve  is  known,  then  (10)  provides  the  means  of  calculating 
the  distribution  of  traction  on  the  crack  faces  that  is  required  to  enforce  the  given  opening. 

The  problem  at  hand,  however,  is  the  inverse  of  this  situation.  That  is,  the  traction 
is  given  along  any  possible  fracture  path  by  the  distribution  of  initial  residual  stress  in 
the  brittle  region  (1).  The  distribution  of  relative  displacement  of  the  crack  faces  along 
the  crack  path  is  unknown  and,  indeed,  the  path  £  itself  is  unknown  a  priori.  Thus,  a 
strategy  is  required  for  determining  the  crack  path  and  the  distribution  of  crack  opening 
required  to  relax  the  initial  residual  stress  in  a  physically  realistic  way.  The  growth  of  the 
crack  negates  the  traction  induced  by  the  initi2J  residual  stress  field  on  the  crack  path. 
This  implies  that  the  crack  opening  is  sought  for  which 


^n(q)  +  <^”(9)  =  0  and  CT,(q)  +  ff®(g)  =  0  along  £  (11) 

In  view  of  (10),  these  conditions  yield  singular  integral  equations  for  the  crack  opening 
distribution  ^i(s)  and  ipy{s). 

A  common  crack  growth  criterion  for  brittle  materials  of  the  kind  imagined  here 
is  that  the  crack  advances  with  Kj  =  Kje  and  Kn  =  0  where  A'/c  the  constant 
fracture  toughness  of  the  material.  Indeed,  the  spiral  shape  of  observed  fracture  paths  in 
consistent  with  this  hypothesis,  and  qualitative  observations  suggest  that  the  condition 
Kji  =  0  may  be  the  precise  mechanism  by  which  the  cracks  search  out  the  spiral  paths 
which  they  follow.  This  idea  is  pursued  in  ib'  \t  section  by  letting  the  crack  seek  a 
path  incrementally  during  growth  according  to  ih;.-  criterion. 


CRACK  PATHS  BASED  ON  K,,  =  0 

The  crack  is  assiimed  to  start  at  the  point  x  =  a,  y  =  0  and  to  grow  in  the  direction 
of  increeising  angle  6.  The  shape  of  the  path  is  described  by  means  of  values  of  the  slope 
of  the  path  0(s)  at  discrete  points.  If  the  arclength  along  the  path  is  measured  from  the 
nucleation  point,  then  the  discrete  p>oints  are  identified  by  the  values  s  =  s*,  fc  =  1, 2, 3, . . . 
8ind  the  values  of  the  slope  at  these  points  are  denoted  by  <t>(3k)  =  It  is  sissumed  that 
the  path  segments  between  those  points  are  circular  arcs,  and  that  the  slope  of  the  path 
is  continuous  at  each  junction  of  adjacent  circular  arcs.  The  slope  of  the  path  between  st 
and  Sk+i  is  lineeir  is  s  and  is  given  by 


i!  \  i  ^ 

<t>(S)  =  <Pk -  +  <Pt+l -  ,  St  <  S  <  S/k+i 

«*+I  St  st+j  St 


(12) 


Fig.  3.  An  increment  of  the  un¬ 
known  crack  path  for  .s*  <  ,s  < 
Hk+t  represented  as  a  circular  arc 
with  known  endpKDints. 


The  radius  of  curvature  p(s)  of  the  path,  constant  within  each  segment,  is 

=  ^'(j)  =  - —  ,  Sk  <  s  <  Sk+i  (13) 

P(«)  ^k+i  -  •»* 

If  the  rectangular  coordinates  of  the  point  located  by  s*  are  (xt,yt),  then  it  is  evident 
from  Figure  3  that  the  coordinates  of  any  point  in  the  interval  st  <  s  <  Si+j  are 

i(s)  =  xt-(-p(s)[sin^(s)-sin«^t]  ,  y(s)  =  yt  -  p(s)  [cos^(s)  -  cos^t]  (14) 

The  quantities  r(s)  and  $(s)  needed  for  the  calculation  are  readily  obtained  from  these 
rectangular  coordinates.  A  similar  description  has  been  used  in  [4]. 

With  r(s),  6(s)  and  ^(s)  determined  explicitly,  the  functions  Eij(9,s)  appearing  in 
(10)  can  be  evaluated  by  following  the  steps  outlined  in  Appendix  II.  These  functions  are 
bounded  and  continuous  for  the  class  of  curves  C  being  considered  here.  With  the  left 
side  of  the  equation  (10)  given  in  terms  of  the  initial  residual  stress  distribution  according 
to  (1),  the  equations  provide  a  pair  of  coupled  singulru'  integral  equations  of  the  first  kind 
for  V’»('S)  and  rl>y(s).  Accurate  approximate  solutions  to  equations  of  this  tyi)e  can  be 
found  by  Gaussian  integration  procedures  [7].  To  this  end,  the  unknown  functions  are 
expressed  in  the  form 


V>i(s)  = 


g»(^) 

^(«+  -■«)(«-«-)  ’ 


9v(s) 

/(S+  -  s){s-3.) 


to  reflect  the  anticipated  square  root  singularities  in  these  functions  at  the  ends  of  the 
crack.  The  functions  5i(s)  and  3,(s)  are  expected  to  be  bounded  and  continuous,  and 
to  be  expressible  to  any  degree  of  accuracy  as  a  finite  series  of  Chebyshev  polynomials  of 


the  first  kind.  The  weight  function  appearing  in  the  orthogonality  relationship  for  these 
polynomials  is 


which  is  provided  by  the  representation  (15).  The  smooth  functions  gz{s)  and  gy{s)  are 
represented  by  a  series  of  the  first  N  Chebyshev  polynomials  with  unknown  coefficients. 
The  equations  (10)  are  then  enforced  at  the  N  —  \  discrete  joints 


Qk  =  5(s+  +s.)-b  |(s+  -s.)< 


ife  =  l . N-l 


and  the  integrals  are  approximated  by  the  appropriate  Gaussian  integration  formulas  [8], 
providing  2N  -  2  equations  for  the  2N  unknown  values  yi(s,)  and  y,(si)  for  »  =  1, . . . ,  N. 


The  Gaussiein  integration  points  are 


Si  =  2(s+  +  3-)  +  §(s+  -  s.)cos  ,  i  =  I,...,  N 

and  the  full  system  of  2N  —  2  equations  is 


(18) 


N 

-T 

N  ^ 
»=1 

N 

-y 

N  ^ 
1=1 


^nx(.9k>^i)9r{^i)  't~  ^ny(^t  i  Sj  )gy(Si ) 

s.  -  qt 

\  S<x(9ti  S|)gx(Si)  4-  ^(y(gt  I  Si)ffy(Si)] 

s,  -  gt 


=  -<T?(9t) 


(19) 


The  remaining  two  equations  are  obtained  from  the  conditions  (8)  that  the  crack  closes 
with  zero  net  dislocation  over  its  entire  length.  Approximation  of  (8)  by  means  of  the 
Gaussian  integration  formula  yields 

N  N 

Y.gx(si)  =  0,  53i7»(s.)=0  (20) 

i=l  i=l 

Although  the  coefficients  of  the  system  of  equations  are  algebraically  complicated,  they 
are  known  explicitly  and  the  equations  are  readily  solved  numerically  by  means  of  the 
Gauss-Jordan  method  or  an  equivalent  method. 

Suppose  that  the  path  has  been  found  for  0  <  s  <  s*.  Then  the  process  of  locating 
p>oint  Sjt+i  by  enforcing  the  condition  Kjj  =  0  is  as  follows.  First,  the  length  of  the 
increment  ds  =  sjt+i  — st  is  specified.  Then,  the  mode  II  stress  intensity  factor  is  calculated 
for  three  trial  selections  of  ^t+i,  namely,  —  e,  =  (t>k  and  + 

e  where  «  is  a  specified  small  number,  typically  about  0.02.  The  three  corresponding 
values  of  mode  II  stress  intensity  factor,  say  A'J";,  Kjj  and  A/},  then  permit  a  parabolic 
interpolation  for  Kji  versus  4>t+\.  The  root  for  Kjj  —  0  provides  afirst  approximation  for 
the  desired  value  of  <t>k+i.  The  value  of  A//  corresponding  to  this  root  is  then  calculated. 
If  the  value  is  less  in  magnitude  than  some  predetermined  magnitude,  typically  0.001  A;, 
then  the  process  is  terminated.  If  this  condition  is  not  fulfilled,  then  the  root  is  used 
as  the  basis  for  a  second  iteration  following  the  template  outlined  above.  The  process  is 
repeated  until  convergence  is  achieved,  whereupon  the  step  to  locate  st+z  is  begun. 

Representative  results  for  the  variation  of  mode  I  stress  intensity  factor  with  arclength 
along  the  computed  crack  paths  eire  shown  in  Figure  4  for  the  case  when  /im  =  ^  and 
k'm  =  in  Figure  5  when  Pm  =  and  Vm  =  v,  and  in  Figure  6  when  fim  =  p/2 
and  Vm  =  V.  The  corresfKjnding  crack  paths  are  also  shown  as  insets  in  the  figures. 
Some  qualitative  features  are  evident  from  the  calculated  shapes  of  the  spiral  cracks.  For 
example,  an  increeise  in  the  stiffness  ratio  causes  the  crack  to  spiral  away  from  the 

inclusion  at  a  faster  rate. 


CONCLUDING  REMARKS 

There  are  no  exact  elasticity  solutions  for  spiral  cracks  in  bi-materitds  systems  of 
the  type  considered  here.  Consequently,  the  accuracy  of  the  numerictd  algorithm  for 
solving  problems  could  be  ascertained  only  by  considering  some  limiting  cases.  When 
the  two  materials  have  the  same  elastic  properties,  that  is  when  a  =  ^  =  0,  there  are 
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Fig.  4.  Mcxle  I  stress  intensity  fac¬ 
tor  veisxis  crack  length  for  the  case 
when  the  two  materials  have  iden¬ 
tical  elastic  properties.  The  inset 
shows  the  computed  crack  path. 


Fig.  5.  Mode  I  stress  intensity  fac¬ 
tor  versus  crack  length  for  fim  = 
2p,  and  the  computed  crack  path. 


Fig.  6.  Mode  I  stress  intensity  fac¬ 
tor  versus  crack  length  for  pm  = 
p/2,  and  the  computed  crack  path. 


several  solutions  available  for  cracks  that  are  straight  or  are  circuW  arcs  concentric  with 
the  interface.  In  such  caises,  the  computed  stress  intensity  factors  agreed  with  the  exact 
results  to  within  0.1%  for  TV  =  15.  There  is  also  an  elastic  crack  solution  for  the  case  of 
distinct  material  properties  but  with  the  crack  lying  in  the  circular  interface  and  extending 
over  an  arbitrary  arc  of  the  interface  [9,10].  However,  the  present  analysis  is  not  applicable 
to  the  case  of  an  interface  crack,  and  this  matter  is  currently  being  pursued. 

The  feature  that  each  curve  for  mode  I  stress  intensity  factor  has  a  local  maximum 
at  a  value  of  s/a  of  about  unity  is  interpreted  in  the  following  way.  If  a  crack  is  nucleated 
at  some  point  on  the  interface  between  the  cylinder  and  the  surrounding  matrix  and  it 


s 


subsequently  grows  along  a  spiral  path  with  more  or  less  constant  fracture  toughness, 
then  that  growth  is  initially  unstable.  For  any  rrark  length  between  r.ero  and  the  length 
when  Ki  is  meiximum,  the  crack  grows  into  a  region  of  increasing  stress  intensity  factor. 
Thus,  growth  is  spontaneous  in  this  region.  Beyond  the  maximum  point  the  growth  is 
stabilized.  The  actual  process  is  probably  more  complex  than  this  simple  idea  suggests, 
perhaps  involving  some  dynamic  crack  growth  phase. 

A  point  of  some  practical  significance  is  that  the  spiral  cracks  continue  to  grow  away 
from  the  inclusion  if  the  inclusion  is  significantly  stiffer  than  the  surrounding  insulator.  On 
the  other  hand,  if  the  inclusion  is  significantly  less  stiff  than  the  insulator,  then  the  radial 
reinge  of  the  spiral  crack  appears  to  be  limited.  Consequently,  a  reduction  of  stiffness  of 
the  inclusion  by  material  selection  or  leaving  a  hollow  core  would  tend  to  diminisli  the 
prospect  of  inclusion-to-inclusion  cracking. 
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APPENDIX  I:  STRESS  FIELD  DUE  TO  A  SINGLE  DISLOCATION  NEAR 
A  CIRCULAR  INCLUSION 

The  configuration  for  which  the  fundamental  dislocation  solution  is  required  is  shown 
in  Figure  7.  The  dislocation  is  located  at  a  point  with  polar  coordinates  ri,6i.  The 
Burgers  vector  of  the  dislocation  heis  recteingular  components  ft,  and  by.  The  mechanical 
field  of  the  dislocation  is  conveniently  described  in  rectangular  coordinates  oriented  so 
that  the  dislocation  lies  on  one  coordinate  axis.  Thus,  the  rotated  ({,t?)  rectanguW 
coordinate  system  is  introduced  as  shown  in  the  figure.  For  any  given  pioint  (x,  y)  in  the 


Fig.  7.  Configuration  which  pro¬ 
vides  the  fundamenta]  elastic  dis¬ 
location  solution  necessary  to  for¬ 
mulate  integrzil  equations  model¬ 
ing  spiral  crack  growth. 


original  coordinate  system  and  any  Burgers  vector  (6*,  6y),  the  corresponding  position  in 
the  rotated  coordinates  is 


{  =  I  cos Si  +  ysinBi,  t]  =  —x sin  6d  +  y  cos 64 
and  the  corresponding  Burgers  vector  is 

b(  =  bi  cos  0i  -f  iy  sin  64,  b,,  =  —bx  sin  9i  -i-  by  cos  6i 


(21) 


(22) 


In  the  rotated  coordinate  system,  the  Airy  stress  function  for  the  stress  distribution 
outside  the  inclusion  has  been  given  by  Dundurs  and  Mura  [11]  as 


^lb^ 


2jr(l 

n{0‘  -  0)  f(Q/T-rf)^(r^/a^  -  l)(2^rrf  -  - 

2(1  +  0) 


x(C,'?)  = 

iTT^r  -  I')  (  1  -  A 

+  )  {(4  -  ra)lnri  +  [(f  -  aVrd)lnr2  -  Onr]  -  -w) 


+ 


2(0-^).  2, 2  111^+^ _ 

1  +  0  ’1-0  (1  -/?)(!  -I- Q- 2/?) 


2ri  J 


where 


r=\/^M^,  n  =  y/a-Tif  +r0  ,  =  ■/({  -  a^/r^)*  -f  r0 

w  =arctan  |r;/^| ,  wi  =  arctan  [q/({  —  r^)] ,  u>2  =  arctan  [t7/(^  —  a^/r^)] 

eind  the  material  parameters  a  and  0  Me  defined  by 


tim(l  -  v)-ti(l  -  I'm) 

Pm(l  •')  +  M(1  I'm) 


0  = 


fi„(l  -  2»/)  -p(l  -  2u„) 
2[lim(l  «')  +  f‘(l  *'m)] 


(23) 


(24) 


(25) 


/O 


a  = 


For  the  special  case  when  the  two  materials  are  the  same,  a  =  /?  =  0  and  the  stress 
function  for  the  region  exterior  t«>  the  inclusion  is 


xU,  n)  =  -  ’’d)  -  ifV]  In  (r?*  +  U  -  rjf)  (26) 

In  either  case,  the  stress  components  referred  to  the  rectangulau'  ^,T7-coordinate  system 
are  derived  from  the  Airy  stress  function  according  to 


-  ^ 

“  dr]^  ~  didr) 


(27) 


Finally,  the  stress  components  in  the  original  rectangular  x,y-coordinate  system  are  ob- 
teiined  by  the  fundamental  transformation  rule  for  the  stress  tensor,  given  in  complex 
notation  by 


2  (^y»  "b  ^ir)  ~  2  >  2(^v»  fi'ir)+*^^iy—  [2(*^')'/  "b  ®  ^  (28) 

The  mathematical  expressions  for  stress  components  as  functions  of  x  and  y  can  be 
obtained  explicitly  for  general  values  of  the  elastic  constants.  However,  the  expressions 
are  too  complex  to  be  either  interesting  or  useful.  Nonetheless,  the  steps  outlined  here 
for  constructing  these  expressions  can  be  implemented  in  a  numericaJ  procedure  in  order 
to  calculate  values  of  these  stress  components  without  approximation. 


APPENDIX  II:  STEPS  IN  CALCULATING  E„x(9,s) 

Consider  a  crack  path  C  in  the  plane  exterior  to  the  circular  inclusion.  The  path  is 
specified  parametrically  by  fioleir  coordinates  r(g),  6{q)  of  each  point  in  terms  of  ^t^clength 
q  along  C  at  the  point.  If  a  dislocation  is  located  at  a  particular  point  defined  by  arclength 
q  =  s  oa  C,  then  the  distribution  of  traction  as  a  function  of  position  q  along  £  is  singular 
(9  —  The  function  E„,(g,  s)  is  the  amplitude  of  this  singularity  for  the  normal 

traction  due  to  the  x-component  of  the  Burgers  vector  of  the  displacement.  The  purpose 
of  this  section  is  to  state  the  steps  involved  in  actually  calculating  this  amplitude  from 
the  dislocation  solution  in  the  previous  section. 

If  the  position  of  the  dislocation  is  specified  by  a  particuleir  value  of  arclength  s  then 
the  parameters  and  of  Appendix  I  are  replaced  by  r(s)  and  0{s),  respectively. 

(i)  The  first  step  in  the  process  is  to  calculate  the  stress  components  <rf((i,r))  ...  for 

arbitrary  Burgers  vector  This  requires  that  the  the  stress  function  (23)  must 

be  differentiated  according  to  (27),  a  task  most  easily  accomplished  by  means  of  a 
computer  program  capable  of  symbolic  mathematical  operations. 

(ii)  To  find  the  particular  component  of  interest  here,  set  6^  =  1  and  by  =  0.  Then  make 
the  substitutions  implied  by  (22)  in  the  expressions  for  stress,  that  is,  substitute 

— »  cos0(s)  and  fr,  — ♦  — sintf(s)  into  (T(((^,t;)  and  the  other  stress  components. 

(iii)  Next,  it  is  necessary  to  treinsform  coordinates  from  the  i,r]  system  to  the  fixed 
coordinates  x,y  in  the  physical  plane.  For  a  particular  field  point  on  £  represented 
by  arclength  q  2ilong  the  curve,  x  =  r{q)cos6{q)  and  y  =  r(g) sin 6(g).  Thus,  the 
coordinate  transformation  indicated  by  (21)  takes  the  form 

i  -4  r(g)  cos[6(g)  -  6(s)] ,  r)  -»  r(g)  sin[6(g)  -  6(s)]  (29) 

(iv)  At  this  point,  the  stress  comjKjnents  Eire  expressed  in  terms  of  the  coordinates  s 
and  g  along  £,  and  the  only  task  remaining  is  to  transform  these  stress  components 


to  the  loccil  iyTi  coordinate  system  on  C  according  to  (4).  For  the  particular  stress 
r<>mp<>iu'nt  of  iuJon'sl,  t\is  is  :irroiijj>lish(*rl  hy 


S„r(iy,s)  = 

(30) 

+  5kr;>/(9.«)  -  «'«(9,s)]cos2[(«6((7)  -  ^(s)]  -  Ci^(q,s)sm2[<{>(q)  -  ^(i)] 

These  steps  ai'e  readily  implemented  in  a  computational  procedure.  It  is  noted  that 
the  implementation  requires  no  approximation  in  the  calculation  of  the  factors  l^,j{q,s) 
for  riny  values  of  q  and  s. 


